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Abstract 



We present a Monte Carlo study of the one-component (f> model on the 
cubic lattice in three dimensions. Leading order scaling corrections are stud- 
ied using the finite size scaling method. We compute the corrections to scaling 
exponent u with high precision. We determine the value of the coupling A at 
which leading order corrections to scaling vanish. Using this result we obtain 
estimates for critical exponents that are more precise than those obtained 
with field theoretic methods. 
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1 Introduction 

The divergence of quantities like the correlation length £ or the magnetic suscepti- 
bility x i n the neighbourhood of a critical point is described by scaling laws 

£ocr" , xoct" 7 , (1) 

where t = \T — T c \/T c gives the distance from the critical point. However, such 
scaling laws are valid in this simple form only in an infinitesimal neighbourhood 
of the critical temperature on infinitely large systems. Monte Carlo simulations 
however are performed with finite systems. Therefore the analysis of the resulting 
data requires knowledge of corrections to scaling. A similar observation holds for 
experimental data of critical systems that are taken at a finite distance from the 
critical temperature. 

While renormalization group Q (see also, e.g., ref. f2|) predicts the structure 
of corrections qualitatively, an understanding on a quantitative level is needed for 
the correct interpretation of Monte Carlo or experimental data. From e-expansion, 
perturbation theory in three dimensions and high temperature series expansions we 
know that leading corrections are proportional to £ _CJ , with u rj 0.8 for the uni- 
versality class of the three dimensional Ising model || [|, |[| |, [7] . The systematical 
error that is quoted for u is about 1% to 5%. In a study || of universal amplitude 
ratios of the Ising universality class the uncertainty of the estimate of u turned out 
to be a major source of systematic errors. Recent Monte Carlo simulations || [L(J 
indicate that the value of u could be considerably larger than 0.8. Little is known 
about sub-leading corrections to scaling. 

Already in refs. || 0] it was suggested that the study of models that interpolate 
between the Gaussian model and the Ising model should be used to study leading 
order corrections to scaling. Such models allow to vary the amplitude of the cor- 
rections to scaling. In particular they allow to eliminate leading order corrections 
by a suitable choice of the parameters of the action (Hamiltonian). This idea was 
recently implemented in the framework of Monte Carlo simulations and finite size 
scaling |I| |ll|. 

In this paper we study the one-component 4 (or Landau-Ginzburg) model in 
three dimensions on a simple cubic lattice. The action is given by 

S = £{-2k£<M>* + a + ^ + A(^-1) 2 } , (2) 

x /a 

where the field variable <p x is a real number and x labels the lattice sites. /i labels the 
directions and \x is an unit- vector in /x-direction. The Boltzmann factor is exp(— S). 
For A = we get the Gaussian model on the lattice. In the limit A = oo the Ising 



model is recovered. Following ref. |Tl| leading order scaling corrections vanish at 



A = 1.0(1). The authors of ref. |0| find A « 1.145. 
The aim of the present study is three-fold: 



We improve the accuracy of the A at which leading order scaling corrections 
vanish. In particular, we give error estimates for the value that is obtained. 

We obtain an accurate estimate of the correction exponent uj. By simulating 
various values of A we are able to vary the strength of the leading order corrections. 

Finally simulations at the optimal A yield accurate results for the critical expo- 
nents v and r). 

In section 2 we discuss how corrections to scaling that arise from the crossover 
from the Gaussian fixed point to the Wilson-Fisher fixed point can be studied by 
finite size scaling. The Monte Carlo algorithm is explained in section 3. In section 
4 we give an overview of the simulations that have been performed. The analysis 
of the data is presented in section 5. In section 6 we compare our results with the 
literature. Finally we give our conclusions and an outlook. 

2 Scaling corrections and finite size scaling 

e-expansion || tells us that leading corrections to scaling are related to the RG-flow 
from the Gaussian fixed point into the Wilson-Fisher fixed point. In Monte Carlo 
simulations of lattices with finite size L the Binder cumulant 

U(L,K,\) = ±% r ± s (3) 

is the most natural quantity to monitor this flow. The magnetization is given by 

m = Y lX 4>x- 

Since we like to study the flow on the critical surface it is useful to consider a 
second phenomenological coupling (i.e. a non-trivial quantity that is invariant under 
RG-transformations). We have chosen the ratio of partition functions with periodic 
and anti-periodic boundary conditions Z a /Z p [|HJ. Note that partition functions 
are by construction conserved under RG-transformations. Therefore also the ratio 
of partition functions Z a /Z p is conserved. 

Instead of computing the Binder cumulant at k c , it is computed at the k for 
that Z a /Zp takes a fixed value on the given lattice. The practical advantage of this 
method is that no errors are introduced by an inaccurate estimate of k c and that 
due to cross-correlations the statistical error of the Binder cumulant at fixed Z a jZ v 
is smaller than that of the Binder cumulant at fixed k. In the following we will 



always fix Z a /Z p = 0.5425, which is according to ref. |10j a good approximation of 



lim Z a /Z p \ Kc . (4) 

L^oo 

The Binder cumulant at a fixed value of Z a /Z p = 0.5425 means: 

U(L,X) = U(L,R(L,\),\) , (5) 

where k(L, A) is determined by 

Z a /Z P (L, R{L, A), A) = 0.5425 . (6) 



The fact that there exists a unique RG-trajectory running from the Gaussian 
fixed point into the Wilson-Fisher fixed point leads to 

U(L,\)=f(a(\)L) + b(\)L-* + ... , (7) 

where we expect x ~ 2. The function a(X) goes to as A — > 0. For some finite 
value of A leading order scaling corrections vanish and a diverges. Reparametrizing 
the scaling function / 

/(c(A) L-») = f(a(X) L) (8) 

and Taylor-expanding yields 

U(L,X) = U* + c 1 (X)L' UJ + c 2 c 1 (X) 2 L-^ + ... . (9) 

In addition to the leading correction L~ u and powers of it we should expect that 
corrections of order L~ 2 , L~ 4 and so on, which exist in the lattice version of the 
Gaussian model, do survive in some form also at the Wilson-Fisher fixed point. 
One example of such a correction is the restoration of the rotational invariance. 
The result of ref. |13j indicates that such corrections exist almost unaltered at the 



Wilson-Fisher fixed point. 

3 The Monte Carlo algorithm 



We followed the idea of Brower and Tamayo [1 1\ and used a combination of the clus- 
ter algorithm and a local Metropolis algorithm for updating the field. We replaced 
the Swendsen-Wang cluster algorithm |15| that was used by Brower and Tamayo 



with the recently proposed wall-cluster algorithm |10[ . The cluster algorithm only 
updates the sign of the field (p. Ergodicity of the update scheme is reached by 
alternating the cluster update with Metropolis updates that allow to change the 
modulus of the field <fi. Below we give the details of the update schemes that were 
used. 

3.1 Metropolis 

In order to make the updating scheme ergodic we performed one Metropolis sweep 
in an update cycle. The proposal for the site x is generated by 

€ = ^ + s(r-^j , (10) 

where r is a random number uniformly distributed in (0, 1] and s parametrizes the 
size of the change. In our study we have chosen s = 3, yielding acceptance rates 
between 0.4 and 0.6, depending on A. 



3.2 Overrelaxation 

In order to speed up the updating of the modulus of the field we added n overre- 
laxation sweeps to the update cycle. A proposal for the field at site x is generated 
by 

<j>' x = 2« J2 cj) y - <j> x , (11) 

y.nn.x 

where y.nn.x means that y is a nearest neighbour of x. This proposal keeps 



1 is ; 
the energy, i.e. d! > with 



-2k^0 x x+a + 



(12) 



constant. This proposal is accepted if the demon variable d G [0, oo) can take over 



d> = d + \(<Pl - If - X^ 2 - If . (13) 

This means that we keep the combined action S com = S + d constant, d is set equal 
to zero at the beginning of the simulation and it is then updated only with the 
overrelaxation updates. 

As a check of the program we measured the expectation value of d: 

/■OO 

< d >= / dxexp(— x) x = 1 . (14) 

Jo 

For example from our simulation with L — 96 at A = 1.1 we obtained < d >= 

1.000003(4). 

The acceptance rate of the overrelaxation step depends on A. For A = the 

acceptance rate is 1 by construction. For A = 1.1, where most of our simulations 

are done, the acceptance rate is still 0.715. 

3.3 Wall-cluster 

We followed the idea of Brower and Tamayo and used the cluster- up date only for 
updating the sign of the field <fi. For that purpose we consider the system as an 
Ising model with a link-dependent coupling constant given by (3 <xy> = 2K\(j) x \\4> y \- 
This results in the link-dependent delete probability 

p d (x,y) = min[l,exp(-2 / 9< a! j / > sign(0 x ) sign(0 y ))] = min[l, exp(-4«0 a; 2/ )] . (15) 

Given the delete probability there is still freedom in the choice of clusters to be 
flipped (i.e. the fields of the cluster are multiplied by —1). In the Swendsen-Wang 
|i~5|| algorithm clusters are flipped with probability 1/2. In Wolff's [16| single cluster 



algorithm the cluster that contains a randomly chosen site of the lattice is flipped. 
Other rules for the selection of the clusters to be flipped can be found in ref. fPT|| . 
In this work we used the wall-cluster algorithm introduced in ref. | i~0| . In this case 



all clusters which intersect with a randomly chosen two-dimensional plane of the 



lattice are flipped. In ref. 10 we show that this version of the cluster algorithm is 



less affected by critical slowing down than the single cluster algorithm. 



4 The Simulations 

We performed simulations at a large range of A values and lattice sizes L. In table |] 
we give an overview of the simulation parameters and the number of measurements 
for each set of simulation parameters. Most of our simulations were performed on 
200MHz Pentium Pro PCs running under Linux. The program is written in C. 
As random number generator we used our own implementation of G05CAF of the 
NAG-library. 

Table 1: Summary of simulation parameters. In the first row we give the value 
of A, in the second row the lattice size L and in the third row the number of 
measurements divided by 3 x 10 6 . 



A 


L 


stat/3 x 10 6 


0.1 


6,8,12,16,24,32,48 


6,5,5,5,4,2,1 


0.2 


6,8,12,16,24,32,48 


6,5,5,5,4,2,1 


0.4 


6,8,10,12,14,16,18,20,24,32 


6,5,5,6,9,6,6,5,5,2 


0.7 


2,3,4,5,6,8,12,16,24 


6,6,6,6,6,6,7,6,7 


0.8 


2,3,4,5,6,8,12,16,24 


6,6,6,6,6,6,5,5,5 


0.9 


2,3,4,5,6,7,8,9,10,12,14,16 


10,10,10,10,30,30,20,20,20,20,20,20 


1.1 


2,3,4,5,6,7,8,9,10,11,12,13,14,15, 


8,8,8,8,25,30,30,32,30,11,10,6,6,6, 




16,18,20,22,24,28,32,40,48,64,96 


5,5,5,5.5,4,6,4,2.5,2.6,1,0.5 


1.145 


2,3,4,5,6,7,8,9,10,11,12,13,14, 


8,8,8,8,25,30,30,30,30,11,10,6,6, 




15,16,18,20,24,32,48,64 


6,5,5,6,4,3,1.7,0.5 


1.3 


2,3,4,5,6,7,8,9,10,12,14,16 


10,10,10,10,30,30,20,20,20,20,20,20 


1.4 


2,3,4,5,6,7,8,12,16,24 


6,6,6,6,6,6,6,7,6,6 


1.5 


2,3,4,5,6,7,8,12,16,24 


6,6,6,6,6,6,5,5,5,5 


2.5 


6,8,12,16,24,32 


6,5,5,5,4,2 



Per measurement we performed one Metropolis sweep, one or two overrelaxation 
sweeps and four or seven wall-cluster updates. The last wall-cluster update in 
this sequence is also used for the measurement of the boundary variable. The 
direction of the plane used for the wall-cluster is taken from a fixed sequence. First 
perpendicular to the 1-direction, then perpendicular to the 2-direction, and then 
perpendicular to the 3-direction and so on. The position of the plane is chosen 
randomly. 

The update of a single site with Metropolis takes 2.0 x 10 _6 sec and the overre- 
laxation 0.4 x 10 _6 sec. The construction of one site of the cluster takes on average 
4.5 x 10~ 6 sec. 

We computed the integrated autocorrelation times of the observables that were 
measured. For example for our largest lattice size L — 96 at A = 1.1 and k ~ 
k c , using one Metropolis sweep, two overrelaxation sweeps and seven wall-cluster 
updates per measurement, we obtain r x = 3.56(2) and r b = 1.93(1) in units of 



measurements for the integrated autocorrelation time of the magnetic susceptibility 
X and the boundary variable, respectively. 

The average sum of the sizes of the clusters per volume that are flipped in 
one update step are fitted with the ansatz S/V = CL X . We obtain S/V = 
1.394(4) L-°- 4874(7) taking into account L = 48, 64 and L = 96 for A = 1.1 and 
K ~ K c . The x 2 /d.o.f. = 4.9 is rather large. Taking only L = 64 and L = 96 yields 
S/V = 1.410(8) L-°- 4901 ( 14 ). 

This result has to be compared with S/V = 1.008(4) L -0 - 527 ^) that was found 



in ref. [[L(| for the standard Ising model on the cubical lattice. It seems that in the 
case of the 4 model there are rather strong corrections to the simple power law. 

The simulations were performed at k s that were the best estimates of R available 
at the start of the simulations. In order to evaluate observables at k values different 
from n s we used second order Taylor expansion in k at k s . We always checked that 
the error made by truncating the Taylor series is well below the statistical error. 
The Taylor coefficients are obtained from the simulations. 

In total we used about 3.5 years of Pentium Pro CPU-time for these simulations. 
The runs for A = 1.1, which turned out to be close to the optimal A, took about 
570 days. 

5 Analysing the data 

5.1 The Binder cumulant and corrections to scaling 

We analysed the Binder cumulant at Z a /Z p = 0.5425 fixed in order to study correc- 
tions to scaling and to find the value of A where leading order corrections to scaling 
vanish. To get a first impression we have plotted our data for A = 0.4, 0.8, 1.1, 1.5 
and 2.5 for 5 < L < 24 in figure 1. For A = 0.4 and A = 0.8, U is decreasing with 
increasing lattice size L. For A = 1.1 starting from about L — 10 the value of U 
stays constant within error-bars. Therefore the value of A at which leading order 
corrections to scaling vanish should be very close to 1.1. Going to A = 1.5 and 
A = 2.5 the value of U is increasing with the lattice size. For a given lattice size U 
is monotonically increasing with decreasing A. 

Next we analyzed our data for U in a more quantitative fashion. In a first 
attempt we fitted the data with the simple ansatz 

U(L,X) =U* + Cl (X) L~ u . (16) 

Our fit results are summarized in table |2|. The data which are included into the fit 
are determined by two criteria. First the lattice size has to be larger than a minimal 
lattice size L > L min . Second the distance d = \U(L, X) — U*\ has to be smaller than 
a maximal distance d < d max . The first criterion allows to control effects of higher 
order corrections in general, while the second specifically controls higher corrections 
that involve the same scaling field as the leading correction (L~ 2u) , L~ 3uJ , ...). 
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Figure 1: The Binder cumulant U at Z a /Z p = 0.5425 as a function of the 
lattice size L for A = 0.4,0.8, 1.1, 1.5 and 2.5. 
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Figure 2: The Binder cumulant U at Z a /Z p = 0.5425 as a function of the 
rescaled lattice size L' = a(X)L for A = 0.1,0.2,0.4,0.7 and 0.8. 



The numbers for x 2 /d.o.f. indicate that at least L m i n = 12 and d max = 0.04 is 
required to obtain a reliable fit. Results from the fits with a x 2 /d.o.f. close to one 
give uj = 0.85 up to 0.87 for the correction to scaling exponent. This value is clearly 
larger than that given by most recent field theoretic work ||. 

Values for ci(A) are not listed here. They will be discussed in more detail below. 

In a second attempt the data were fitted with the ansatz 



U(L,\) = U* + Cl (\) L— + c 2Cl (A) 2 L 



-2w 



(17) 



The additional term requires only one further parameter in the fit. Our results are 
summarized in table ||. We note that, as expected, data with much larger values of 
d = U — U* can be included in the fit. The variation of uj with varying L min and 
d m ax is considerably reduced compared with the previous fit. 

As our final result we quote uj = 0.845(10), where we take into account the 
variation of the result with L min and d max . 

For a small number of fits we give in table [| the results for Ci(A). We see that 
for A = 1.1, within error-bars, the leading corrections to scaling vanish. Converting 
the error-bar of Ci(A) to A we get X opt = 1.100(7) from L min = 12, X opt = 1.102(8) 
from L min = 14, and X opt = 1.095(12) from L min = 16. X opt is the A where leading 
order corrections vanish. 



Table 2: Fit results for the Binder cumulant evaluated at Z a /Z p = 0.5425 fixed. 
The ansatz is given in eq. (|7^). We give results for various minimal lattice sizes 
L m i n and maximal distances d max of U from its fixed point value U*. uj is the 
correction to scaling exponent. 



J-'min 


U"rnax 


rVd.cf. 


U* 


UJ 


8 


0.06 


5.36 


1.60351(8) 


0.865(3) 


8 


0.04 


2.45 


1.60300(8) 


0.866(5) 


8 


0.03 


2.06 


1.60286(8) 


0.866(7) 


12 


0.06 


2.02 


1.60398(13) 


0.840(5) 


12 


0.04 


1.43 


1.60370(14) 


0.856(6) 


12 


0.03 


1.21 


1.60344(15) 


0.864(10) 


16 


0.06 


1.85 


1.60435(20) 


0.817(8) 


16 


0.04 


1.31 


1.60379(22) 


0.853(11) 


16 


0.03 


1.18 


1.60362(23) 


0.863(12) 



In order to check the existence of the scaling function of eq. (0) we have plotted 
in figure 2 U for A = 0.1, 0.2, 0.4, 0.7, and 0.8 as a function of 



U 



' c(X) 
,c(0.1) 



-l/w 



Table 3: Fit results for the Binder cumulant at Z a /Z p = 0.5425 fixed. Compared 
with the previous table a term proportional to L~ 2uJ has been included into the 
ansatz (|T7p. 



l-'rniri 


Umax 


X 2 / d -cf- 


U* 


UJ 


c 2 


8 


all 


3.08 


1.60318(8) 


0.847(2) 


-1.01(2) 


8 


0.07 


1.83 


1.60276(9) 


0.847(3) 


-1.64(7) 


8 


0.05 


1.94 


1.60277(9) 


0.849(4) 


-1.57(15) 


12 


all 


1.19 


1.60355(13) 


0.839(5) 


-1.17(5) 


12 


0.07 


0.95 


1.60328(15) 


0.846(6) 


-1.57(10) 


12 


0.05 


1.01 


1.60328(16) 


0.846(6) 


-1.54(29) 


16 


all 


1.02 


1.60363(21) 


0.839(9) 


-1.31(10) 


16 


0.07 


0.95 


1.60344(23) 


0.845(11) 


-1.65(18) 


16 


0.05 


1.02 


1.60336(26) 


0.845(11) 


-1.97(48) 



The values for c are taken from table [|. We see that almost all data-points fall 
nicely on a unique curve. Only for the smallest lattice sizes that have been included 
(L = 6) a small deviation is visible. 

Finally we studied differences of the Binder cumulant at Z a /Z p = 0.5425 fixed 
for different values of A. We define 



AU(L,\ 1 ,X 2 ) = U(L,\ 1 )-U(L,\ 2 ) . 



This way U* is cancelled and we expect 



Atf(L,Ai,A 2 ) = [ci(Ai)-ci(A 2 )] L~" + 



(19) 



(20) 



As an example we give in table |5] the results of fits for Ai = 0.9 and A 2 = 1.3. The 
interesting fact is that already starting from L min = 4 the data are well fitted by 
the simple ansatz. Also the value for uj obtained this way is compatible with the 
result obtained above. This means that sub-leading corrections depend very little 
on A and are, to a large extend, cancelled in AU . Hence there should be a good 
chance to study these corrections using perturbation theory. 

Here we try to obtain a better understanding of the sub-leading corrections by 
looking at the results for small L at A = 1.1 in more detail. Corrections to scaling 
vanish very rapidly. Starting from L = 9 the result for the Binder cumulant at 
Z a /Z p fixed is consistent with U* within the error-bars. Fitting 



U(L) -U* = c L- x 



(21) 



we find x ~ 6 ± 1 when including lattices of size L = 4 up to L = 8. It is quite 
surprising that the numerically dominant corrections are governed by such a large 
exponent and not, as one might expect, by x ~ 2. 
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Table 4: The correction to scaling amplitude ci(A) as a function of A from fits 
with the ansatz (|T7f) . We give tie results for three values of L min = 12, 14 and 16. 



0.07 in all three cases. 



A 


J-'min J-^ 


■L-'min J-"^ 


■L-'min -J-0 


0.1 


1.452(44) 


1.454(61) 


1.465(81) 


0.2 


0.809(21) 


0.810(28) 


0.812(37) 


0.4 


0.3882(85) 


0.387(12) 


0.387(16) 


0.7 


0.1485(39) 


0.1482(52) 


0.1462(72) 


0.8 


0.0985(32) 


0.0982(42) 


0.0988(59) 


0.9 


0.0628(22) 


0.0623(31) 


0.0602(41) 


1.1 


0.0001(18) 


0.0006(25) 


-0.0015(34) 


1.145 


-0.0116(16) 


-0.0128(22) 


-0.0139(31) 


1.3 


-0.0479(13) 


-0.0481(16) 


-0.0482(23) 


1.4 


-0.0684(18) 


-0.0686(19) 


-0.0739(30) 


1.5 


-0.0857(20) 


-0.0859(22) 


-0.0896(32) 


2.5 


-0.1933(28) 


-0.1934(36) 


-0.1944(49) 


oo 


-0.3112(39) 


-0.3113(58) 


-0.3125(76) 



Also the authors of ref. |18J found that a correction term with a rather large 
exponent has to be included into the fit ansatz in order to fit the Binder cumulant 
on small lattices. They suggest that this correction exponent is given by 2y^, 
where y^ ~ 2.48 is the RG-exponent related with the external field. To check this 
hypothesis we also looked at the corrections of Z a /Z p at k c . It turns out that they 
vanish as rapidly as the corrections of the Binder cumulant. Hence the observed 
corrections are related to an irrelevant RG-exponent. 



Table 5: Fits of the difference of the Binder cumulant at Ai 
with the ansatz (|2C[). We give dc = ci(Xi 
exponent uo. 



0.9 and A 2 = 1.3 
Ci(A 2 ) and the correction to scaling 



J-'min 


xVd.o.f. 


dc 


us 


2 


2.84 


0.1149(7) 


0.870(4) 


3 


2.83 


0.1132(12) 


0.862(6) 


4 


1.18 


0.1078(18) 


0.838(8) 


5 


1.26 


0.1065(23) 


0.833(11) 


6 


1.39 


0.1076(28) 


0.837(12) 
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5.2 The critical line k c (\) 

In order to obtain a result for the critical coupling k c we fitted our data for the 
ratio of partition functions Z a /Z p with the ansatz 

Z a /Z P (L, k, A) = ZJZ; + 5f|Zf±£ (« _ „ c ) . (22) 

For A = 1.1 we obtain with L min = 12 the result Z a /Z* = 0.54243(9) and 2k c = 
0.3750965(4). Increasing the minimal lattice size that is included in the fit to 
L min = 16 we obtain Z a /Z* p = 0.54244(14) and 2k c = 0.3750966(4). In both cases 
X 2 /d.o.f. is a little larger than 1/2. 

The result Z a /Z* = 0.54244(14) is consistent with previous estimates Z a /Z* = 
0.54334(26){27} from the standard Ising model and Z a /Z* = 0.54254(14) [6]{14} 
from the spin-1 Ising model [Kj. The numbers given in the second and third 



bracket are estimates of systematic errors caused by subleading corrections. In this 
work we skipped a detailed analysis of systematic errors. However it is reasonable 
to assume that systematic errors are of similar size as for the spin-1 Ising model. 
The estimate 2k c = 0.3750966(4) for A = 1.1 will be useful in forthcoming studies 
of scaling laws. 

For the other values of A we determined n r from 



><c 



8ZJZ P 
8k 



ZJZ P + — ^—^ (k c -k) = 0.5425 . (23) 



The results are 

[0.1,0.3734095(13)], [0.2,0.3884251(14)], [0.4,0.3975837(13)], [0.7,0.3925302(20)], 
[0.8,0.3887757(21)], [0.9,0.3845113(36)], [1.145,0.3728926(7)], [1.3,0.3652233(33)], 
[1.4, 0.3602789(19)], [1.5, 0.3553854(21)], and [2.5, 0.3134347(16)] for the pairs [A, 2« c 

5.3 The magnetic susceptibility 

As definition of the magnetic susceptibility we used 




X=T* £<M ■ (24) 



Following refs. ]TD, [L9| we tried to extract the exponent rj from the magnetic sus- 
ceptibility evaluated at Z a /Z p = 0.5425 or at U — 1.6034 fixed. In the following 
the magnetic susceptibility evaluated at Z a /Z p or at U fixed is denoted by % We 
fitted our data with the ansatz 

x(L)=c + dL 2 ^ , (25) 

where c is an analytic correction. Note that also corrections that decay like L~ x 
with x ~ 2 are effectively parametrized by this ansatz. 
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Table 6: Results for rj from fitting the magnetic susceptibility at Z a /Z p 
fixed with the ansatz (\2t 



0.5425 



J-'min 


X 2 /d.o.f. 


V 


d 


c 


4 


25.05 


0.0392(1) 


0.9856(3) 


-0.724(4) 


6 


2.73 


0.0371(1) 


0.9789(4) 


-0.577(8) 


8 


1.23 


0.0364(2) 


0.9765(7) 


-0.504(18) 


10 


1.23 


0.0361(3) 


0.9755(9) 


-0.46(3) 


12 


1.33 


0.0358(4) 


0.9745(14) 


-0.40(7) 


14 


1.29 


0.0357(5) 


0.9740(18) 


-0.35(11) 


16 


1.00 


0.0358(6) 


0.9742(21) 


-0.35(17) 



Table 7: Results for rj from fitting the magnetic susceptibility at U = 1.6034 hxed 



with the ansatz (25|). 



l-'rnin 


xVd-o.f. 


V 


d 


c 


4 


3.24 


0.0345(2) 


0.9701(5) 


-0.31(1) 


6 


1.24 


0.0356(3) 


0.9736(7) 


-0.38(1) 


8 


1.34 


0.0359(4) 


0.9744(12) 


-0.41(3) 


10 


1.28 


0.0365(5) 


0.9767(17) 


-0.50(6) 


12 


1.09 


0.0370(7) 


0.9784(25) 


-0.59(12) 


14 


1.02 


0.0368(9) 


0.9768(32) 


-0.45(20) 


16 


1.26 


0.0366(10) 


0.9768(39) 


-0.45(31) 



First we performed fits for A = 1.1, where leading corrections to scaling vanish. 
Results for various L min are given in table ^| for Z a /Z p = 0.5425 fixed and table [F] 
for U = 1.6034 fixed. As before, L m i n gives the smallest lattice size that has been 
taken into account for the fit. The x 2 /d.o.i. becomes about 1 for L m i n = 8 for 
Z a /Z p fixed, and for L min = 6 for U fixed. 

As a check we performed fits without analytic part c for Z a /Z p = 0.5425 fixed 



x{L) = dL 



2-r, 



(26) 



The results are summarized in table |^. Starting from L min = 28 the result for n 
becomes consistent with the results obtained from the fit eq. (55) that includes 
analytic corrections. 

Finally we checked for systematic errors due to residual leading order corrections 
to scaling at A = 1.1. For that purpose we fitted the magnetic susceptibility at 
Z a /Zp = 0.5425 for A = 0.4,0.8, 1.5, and 2.5 with L min = 6 using the ansatz (p5|). 
The results are n = 0.0493(4), 0.0407(6), 0.0353(5), and 0.0330(5) for the four values 
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Table 8: Results for rj from fitting the magnetic susceptibility at Z a /Z p = 0.5425 

fixed with the ansatz 

corrections. 



Note that the ansatz (\2q) contains no term for analytic 



L-'min 


X 2 /d.o.f. 


V 


d 


12 


3.92 


0.03357(15) 


0.9670(4) 


14 


2.10 


0.03431(18) 


0.9688(5) 


16 


1.35 


0.03467(22) 


0.9700(7) 


20 


1.30 


0.03497(28) 


0.9710(9) 


24 


1.32 


0.03534(37) 


0.9724(13) 


28 


1.60 


0.03545(43) 


0.9728(15) 



of A, respectively. From the results for A = 0.8 and A = 1.5 we get the estimate 



An, 



'eff 



AA 



-0.01 



(27) 



where ry e // denotes the numerical result for r\ obtained from fitting data for lattice 
sizes L = 6 up to L = 24. 

From the previous section we know that the difference of 1.1 and the values 
of A where leading order corrections vanish exactly should be smaller than 0.02. 
Therefore the systematic error in our final estimate of n due to residual leading 
order corrections should be smaller than 0.01 x 0.02 = 0.0002. Note that the lattice 
sizes used to obtain our final result range from L = 8 up to L = 96. 

As final estimate for n we take the result from fitting the magnetic susceptibility 
at Z a /Zp = 0.5425 with the ansatz (|25| ) and L m i n = 12 



n = 0.0358(4) [5] . 



(28) 



The estimate of the systematic error is given in the second bracket. It is obtained 
from the comparison with the ansatz ( p6|) without analytic corrections and from 
the discussion above on L~ w corrections. 



5.4 The exponent 



v 



We compute the exponent v from the slope of the Binder cumulant and the slope 

of Z a /Zp at Z a /Zp = 0.5425 fixed. We fitted our data for A = 1.1 with the simple 

power law ansatz 

~7>TT 

c L l ' v . 



dU 



(29) 



In table |9] we give our fit results for various values of the minimal lattice size L min 
that is included in the fit. The xVd.o.f. becomes reasonably small for L min = 12. 
Also the result for v stays stable when L min is further increased. 
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In table [TO] we give the analogous results for the slope of Z a /Z p . We see that 
the x 2 /d.o.f. are much larger than for the slope of the Binder cumulant. The value 
obtained for v is increasing with increasing L m i n . 



Table 9: Fits of the derivative eq. of Binder cumulant U at Z a /Z p = 0.5425 fixed. 
The ansatz ( |2^ ) is used. As results we give the exponent of the correlation length 
v and the constant c. 



J-'min 


X 2 /d.o.f. 


V 


2c 


6 


2.10 


0.6289(1) 


-1.744(1) 


8 


1.70 


0.6292(1) 


-1.748(2) 


10 


1.76 


0.6293(2) 


-1.750(2) 


12 


1.48 


0.6296(3) 


-1.753(3) 


14 


1.50 


0.6299(3) 


-1.757(4) 


16 


1.76 


0.6299(4) 


-1.756(6) 



Table 10: Fits of the derivative of Z a /Z p at Z a /Z p = 0.5425 fixed. The ansatz 
(p%[) is used. As results we give the exponent of the correlation length v and the 
constant c. 



J-'min 


rVd.o.f. 


V 


2c 


6 


117.0 


0.6254(1) 


-1.1524(4) 


8 


15.1 


0.6273(1) 


-1.1673(5) 


10 


4.9 


0.6281(1) 


-1.1742(7) 


12 


2.9 


0.6286(1) 


-1.1794(12) 


14 


2.2 


0.6290(2) 


-1.1830(15) 


16 


1.9 


0.6292(2) 


-1.1858(20) 



Also for the spin 1 Ising model the authors of ref. [|IIJ found that the slope of 
the Binder cumulant is scaling much better than the slope of Z a /Z p . 

Therefore we take the result obtained from fitting the slope of the Binder cu- 
mulant with L min =12 v = 0.6296(3) as our final result. In order to check for 
systematic errors due to residual L~ w corrections we fitted our data for the slope 
of the Binder cumulant at A = 0.4,0.8,1.5, and 2.5 with L min = 8. We obtain 
v = 0.6363(3), 0.6325(5), 0.6271(4), and 0.6241(4) for the four values of A, respec- 
tively. From A = 0.8 and 1.5 we obtain 



Ai/, 



eff 



AA 



-0.01 



u e ff means here the exponent v obtained from fitting lattices of size L 
L = 24 with the simple power law eq. (|29|) . 



(30) 



8 up to 



15 



From the previous section we know that the difference of 1.1 and the values 
of A where leading order corrections vanish exactly should be smaller than 0.02. 
Therefore the systematic error in our final estimate of v due to residual leading 
order corrections should be smaller than 0.01 x 0.02 = 0.0002. Also note that the 
lattice sizes used to obtain our final result range from L = 12 up to L = 96. We 
arrive at the final result 

z/ = 0.6296(3)[4] , (31) 

where the second bracket gives an estimate of systematic errors. It is obtained from 
the discussion on residual L~ w corrections and from the comparison of fits of the 
slope of the Binder cumulant and of the slope of Z a /Z p . 

6 Comparison with the literature 

There exist exhaustive compilations of results for critical exponents in the literature. 
See for example refs. [^ |18|. In table [Tl] we give only most recent results that reflect 
the state of the art. In ref. [|H|] the spin-1 model was simulated at parameters where 
leading order corrections vanish (up to numerical uncertainties). The approach to 
compute critical exponents is very similar to the one of the present work. The 
results for v and rj are consistent with the results of this study. The error-bars are 



of similar size. In ref. [[T(| the authors found indications that the value of u should 
be larger than that obtained by field-theoretic methods, however they were not able 
to give reliable error estimates. Simulating the (f) 4 model at various values of the 
coupling constant A allowed us to vary the strength of the corrections to scaling. 
Therefore we were able to give a reliable error-estimate for to. 

In ref. f| the standard Ising model was simulated. The authors also use finite 
size scaling techniques to compute the exponents. In the analysis of the data leading 
order corrections to scaling are taken into account. The final result for v is perfectly 
consistent with our present result. The error-bars are a little larger than ours. 
The value for 7] is larger than ours but still compatible when the statistical and 
systematical errors are taken into account. Their value for the correction exponent 
uj = 0.87(9) is consistent with ours but the error-bar is nine times larger than ours. 

Refs. P, ITDJ as well as the present work employ finite size scaling methods that 



were pioneered by Binder [20|. It seems to us that such an approach is more robust 



than the so called Monte Carlo renormalization group method [21]. 

The e-expansion was invented by Wilson and Fisher ||. Most recent results 
obtained from the e-expansion and from perturbation theory in 3 dimensions are 
given in ref. ||. The results for v and i] obtained from the e-expansion are in 
perfect agreement with our results. However the error-bars are considerably larger 
than those obtained now from Monte Carlo. The result for uj is smaller than ours 
but still consistent within error-bars. Parisi 0] proposed to perform perturbative 
expansions directly in three dimensions (3D PT). The result given in ref. |J for v 
is consistent with ours. The value of rj is smaller than ours, but still within the 
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quoted error-bars. In both cases our error-bars are considerably smaller. The value 
of the correction exponent uo obtained from 3D PT || is considerably smaller than 
our result. The results are inconsistent when the error-bars that are quoted are 
taken into account. 

The analysis of high temperature expansions was the first theoretical method 
that produced reliable, non-mean-field values for critical exponents of three dimen- 
sional systems. For a review of early work see ref. |J22| . Still this method gives 



results that are competitive in accuracy with Monte Carlo and field theoretic meth- 



ods. The authors of ref. |23| claim to have the most accurate results obtained from 
high temperature series expansions of the Ising model on the bcc lattice. Their 
result for v is larger than ours. The quoted error-bars do not overlap by a small 
margin. The authors do not quote a value for r\ in their paper. For the convenience 
of the reader we converted 7 = 1.2384(6) given for the exponent of the magnetic 
susceptibility using the scaling relation 77 = 2 — 7/V. (Converting our numbers we 
get 7 = 1.2367(11).) One should note that in this work u> as well as (3 C is taken as 
external input for the analysis. 

We also give results obtained from the series expansion of two parameter models 
that interpolate between the Gaussian and the Ising model, similar to the model 
considered in the present work. While the results of Nickel and Rehr |7j are in 
perfect agreement with our numbers (with larger error-bars) the value obtained by 
Chen, Fisher and Nickel || for v is clearly larger than ours. Also the value 7 = 
1.237(2) obtained by Nickel and Rehr is consistent with ours, while 7 = 1.2395(4) 
of Chen, Fisher and Nickel seems too large compared with our result. 

The values obtained by Nickel and Rehr as well as by Chen, Fisher and Nickel 
for the correction to scaling exponent uo are consistent with our result. However the 
error-bars are considerably larger than ours. 

Further methods to compute critical exponents ("exact renormalization group", 
"coherent anomaly method", ...) are discussed in the literature. These methods 
tend to give less accurate results and are therefore not discussed here. 

Experimental results have been obtained for binary mixtures of fluids, vapor- 
fluid systems and uni-axial (anti-)ferromagnetic systems. In general the results 
agree with the predictions yielded by the theoretical methods discussed above. How- 
ever the results for the exponents are less accurate than the theoretical ones. To 
give only a few examples: The study of the heat capacity of a aniline-cyclohexane 
mixture gives a = 0.104(11) J23. Using the scaling relation dv = 2 — a we get 
v = 0.6320(37). A neutron scattering measurement of an antiferromagnetic FeF 2 
system gave v = 0.64(1) and 7 = 1.25(2) p5|. For a large collection of experimental 
results see ref. ill 
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Table 11: Recent results for critical exponents obtained with Monte Carlo simu- 
lations (MC), e-expansion, Perturbation-Theory in three dimensions (3D,PT) and 
High temperature series expansions. When only v and 7 are given in the refer- 
ence we computed n with the scaling law. These cases are indicated by * . For a 
discussion see the text. 



Ref. 


Method 


V 


77 


UJ 


present work 


MC 


0.6296(3)(4) 


0.0358(4)(5) 


0.845(10) 


m 


MC 


0.6298(5) 


0.0366(8) 




E 


MC 


0.6294(5)(5) 


0.0374(6)(6) 


0.87(9) 


§ 


3D,PT 


0.6304(13) 


0.0335(25) 


0.799(11) 


§ 


e,bc 


0.6305(25) 


0.0365(50) 


0.814(18) 


B 


e,free 


0.6290(25) 


0.0360(50) 


0.814(18) 


23| 


HT 


0.6308(5) 


0.0368(18)* 




1 


HT 


0.632(1) 


0.0388(32)* 


0.854(80)* 


§ 


HT 


0.6300(15) 


0.0365(56)* 


0.825(50)* 



7 Conclusion and outlook 

In this study we determined with high precision the value of the coupling A for 
which leading corrections to scaling vanish. This allowed us to obtain results for 
critical exponents with an accuracy better than that of field theoretic methods. 

The knowledge of the optimal A can also be used for the study of other universal 
properties of the model, e.g. universal amplitude ratios, the effective potential. 

The programme to eliminate leading order corrections can also be applied to the 
4 model with more than one component of the field. The case of two components 
is of particular interest since it is supposed to be in the same universality class as 
superfluid helium systems. For this system there exist experimental results for v 
|26[| that are by far more accurate than the existing theoretical predictions. 

It would be desirable to extend the programme to subleading corrections. How- 
ever these corrections are by far less well understood than leading corrections to 
scaling. Therefore it is unclear how many and what kind of terms should be included 
into the action. 

However, even if it is not possible to remove higher order corrections to scaling 
completely it would be a good check of the reliability of the results to redo the study 
with a different lattice or with more than a next to nearest neighbour coupling. 



8 Acknowledgements 

We thank K. Pinn and T. Torok for discussions and comments on the manuscript. 



18 



References 



[i] 

[2] 

[3] 
[4] 

[5] 
[6] 

[7] 
[8] 

[9] 



F.J. Wegner, Phys.Rev.B 5 (1972) 4529. 

J. Cardy, Scaling and Renormalization in Statistical Physics, Cambridge Uni- 
versity Press 1996. 

K.G. Wilson and M.E. Fisher, Phys.Rev.Lett. 28 (1972) 240. 

G. Parisi, J.Stat. Phys. 23 (1980) 49. 

R. Guida and J. Zinn- Justin, |cond-mat/980324q , J.Phys.A 31 (1998) 8103. 

J.H. Chen, M.E. Fisher and B.G. Nickel, Phys.Rev.Lett 48 (1982) 630; 
M.E. Fisher and J.H. Chen, J.Physique (Paris) 46 (1985) 1645. 

B.G. Nickel and J.J. Rehr, J. Stat. Phys. 61 (1990) 1. 

M. Caselle and M. Hasenbusch, J.Phys.A 30 (1997) 4963. 

H. G. Ballesteros, L. A. Fernandez, V. Martin-Mayor, G. Parisi, and J. J. 
Ruiz-Lorenzo, |cond-mat/9805125| , J. Phys. A 32 (1999) 1. 



[10] M. Hasenbusch, K. Pinn and S. Vinti, 

Critical Exponents of the 3D Ising Universality Class From Finite Size Scaling 
With Standard and Improved Actions, 



hep-lat/9806"0i~2l to be published in Phys.Rev.B . 



[11] H.G. Ballesteros, L.A. Fernandez, V. Martin-Mayor and A. Munoz-Sudupe, 
|hep-lat/ 9805021 , Phys.Lett.B 441 (1998) 330. 



[12 
[13 

[14 
[15 
[16 
[17 
[18 
[19 



M. Hasenbusch, Physica A 197 (1993) 423. 

M. Campostrini, A. Pelissetto, P. Rossi, and E. Vicari, |cond-mat/970508€ 



Phys.Rev.E 57 (1998) 184. 

R.C. Brower and P. Tamayo, Phys.Rev.Lett. 62 (1989) 1087. 

R.H. Swendsen and J.-S. Wang, Phys.Rev.Lett. 58 (1987) 86. 

U. Wolff, Phys.Rev.Lett. 62 (1989) 361. 

W. Kerler, Phys.Rev.D 47 (1993) 1285. 

H.W.J. Bldte, E. Luijten and J.R. Heringa, J. Phys. A (1995) 6289. 

H. G. Ballesteros, L. A. Fernandez, V. Martin-Mayor, A. Munoz-Sudupe, G. 
Parisi, and J.J. Ruiz-Lorenzo, Phys.Lett. B 400 (1997) 346. 

19 



[20] K. Binder, Z.Phys. B 43 (1981) 119; K. Binder, Phys.Rev.Lett. 47 (1981) 693. 

[21] S.K. Ma, Phys.Rev.Lett. 37 (1976) 471; G.S. Pawley, R.H. Swendsen, D.J. 
Wallace and K.G. Wilson, Phys.Rev. B 29 (1984) 4030. 

[22] M.E. Fisher, Rep. Prog. Phys. 30 (1967) 615. 

[23] P. Butera and M. Comi, |hep-lat/9703"0l8| , Phys.Rev. B 56 (1997) 8212. 

[24] PF. Rebillot and D.T. Jacobs, J.Chem.Phys. 109 (1998) 4009. 

[25] D.P Belanger and H. Yoshizawa Phys.Rev.B 35 (1987) 4823. 

[26] J. A. Lipa, D.R. Swanson, J. Nissen, T.C.P. Chui and U.E. Israelson, 
Phys.Rev.Lett. 76 (1996) 944; 
L.S. Goldner, N. Mulders and G. Ahlers, J. Low Temp.Phys. 93 (1992) 131. 



20 



